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Abstract 

We construct the hydrodynamic equations for suspensions of self-propelled particles (SPPs) 
with spontaneous orientational order, and make a number of striking, testable predictions: (i) 
SPP suspensions with the symmetry of a true nematic are always absolutely unstable at long 
wavelengths, (ii) SPP suspensions with polar, i.e., head-tail asymmetric, order support novel 
propagating modes at long wavelengths, coupling orientation, flow, and concentration, (iii) In a 
wavenumber regime accessible only in low Reynolds number systems such as bacteria, polar-ordered 
suspensions are invariably convectively unstable, (iv) The variance in the number N of particles, 
divided by the mean (N), diverges as (iV) 2//3 in polar-ordered SPP suspensions. 

PACS numbers: 89.20.-a, 05.40.-a, 47.90. +a, 47.35. +i 
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Fish, birds, and swimming cells [l], |2|, [3| are frequently found to move coherently in large 
groups ||, ||, || 0, H through the fluid medium they inhabit. Can general principles, like 
those so used so successfully for ordered phases at equilibrium || , reveal the laws governing 
the long-wavelength dynamics and fluctuations of these striking and ubiquitous examples 
of liquid-crystalline |TD[ order in suspensions of self-propelled particles (SPPs)? Although 
we know of no physics experiments on ordered SPP suspensions, Gruler jTIJ has studied 



ordered phases of living cells on a solid substrate, and finds "living liquid-crystalline" fl2 



phases corresponding to two distinct types of cells: apolar, that is, elongated but head-tail 
symmetric, and polar, distinguishing front from rear. The ordered phases formed by polar 
SPPs have a nonzero macroscopic drift velocity Vq. Those formed by apolar SPPs have a 



macroscopic axis n of orientation but, like true nematics [13], do not distinguish n from 
— n; they thus have vq = 0. Migratory cells such as white blood cells are polar in this 
sense, while melanocytes, which distribute pigment in the skin, are apolar pT| . Although 
an isolated melanocyte is incapable of directed motion, we still term it an SPP because 
it displays spontaneous, self-generated energy-dissipating activity, in the form of a pair of 
symmetrically pulsating dendrites [|ll|]. The pioneering studies of [§] and || |7], [| consider 
ordered flocks drifting through a passive frictional background. This is appropriate for polar 
SPPs on a substrate as in [ll], but not for SPPs immersed in a bulk fluid. In this Letter 
we use symmetry and conservation laws to construct the complete equations of motion for 
small, long-wavelength disturbances in polar as well as apolar ordered SPP suspensions, 
including the flow of the ambient fluid. 

Our main results, which we now summarize, contain striking, experimentally testable 
signatures of the nonequilibrium nature of ordered SPP suspensions, and of the crucial 
role of the hydrodynamics of the ambient fluid medium: Purely nematic order in SPP 
suspensions is always destabilized at small enough wavenumber q, by a coupled splay of the 
axis of orientation and a corresponding Taylor-Couette-like circulation of the velocity field, 



oriented near 45° to the nematic axis, with a growth rate linear in q [14 



Polar ordered suspensions display some rather original propagating modes as a result of the 
interplay of hydrodynamic flow with fluctuations in the ordering direction and the concen- 
tration: a pair of bend-twist waves (with no analogue in the work of |J), with wavespeeds as 
in see Fig. [I] and eqn. (|7|), and three waves which are a combination of splay, concentration 
and drift (a generalization of those in II), whose speeds are better understood from Fig. 
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FIG. 1: A polar plot of speed of the bend- twist modes as a function of the angle 9 between the 
propagation direction and the z axis. 

than from an equation. 

The results in Figs, p] and ||] and Eqn. (0) are obtained when viscous damping is ne- 
glected, which is valid for wavenumbers q <C vq/u where v is a typical kinematic viscosity. 
Experiments on bacterial suspensions are, however, likely to be in the Re <C qa <C 1 regime, 
where Re = v^ajv is the Reynolds number of an SPP of size a. In that (Stokesian) limit 
we find, remarkably, that a polar-ordered suspension is always unstable for wavevectors q 
near 45° to the nematic axis, with a growth rate ~ v^/a, independent of |q|. The instability 
is "convective" : it travels with a speed ~ vq as it grows. Lastly, number fluctuations in 
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Sound speed ( in units of v Q ) 



FIG. 2: A polar plot of the speed of the splay-concentration and drift waves as a function of the 
angle 6 between the propagation direction and the z axis. 

polar-ordered SPP suspensions are anomalously large. The variance ((SN) 2 ) in the number 
of particles, scaled by the mean N, is predicted to diverge as N 2 ^ 3 

The results above follow from a set of hydrodynamic equations of motion, which we now 
construct by generalizing |J. Assume we have nematic or polar ordered phases of SPPs, with 
orientation given by a unit director field n aligned on average in the z direction. The nematic 
is invariant under n — > — n; the polar-ordered state is not. The slow variables M are (a) the 
"conserved modes" which, for an incompressible suspension are |15j the fluctuations Sc(r, t) 
at point r and time t in the local concentration c of suspended particles about its mean c 
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and the total (solute + solvent) momentum density g(r, t) = pu(r, t), where p is the constant 
mass density of the suspension and u the hydrodynamic velocity field; and (b) the "broken- 
symmetry" variables which, for both polar and nematic phases, are simply the deviations 
<5n_|_ = n — z. In polar-ordered suspensions, the drift velocity of SPPs relative to the fluid 
is t>on. Since these are driven systems, the forces entering the equations of motion need 
not arise from a free-energy functional, and any term not explicitly forbidden by symmetry 
or conservation laws must be included. Moreover, there are a priori no relations amongst 
the phenomenological parameters other than those required by the geometrical symmetries 
of the problem. Note that our coarse-grained description applies on length-scales ^> a, the 
particle size. 

We now present the equations of motion which follow from these arguments, and explain 
the physical origin of the terms therein. We begin with the polar-ordered case. The equations 
for the nematic case follow by dropping terms which violate z — *> —z symmetry. The director 
field obeys 

d t 5nj_ = — Aii;o<9 2 <5n_i_ — o~\ V±5c + 

^(d z u_L - V±u z ) + -7 2 (<9*ui + V±_u z ) + ... (1) 

The first term on the right-hand side of (|l|) represents advection of distortions by the mean 
drift v ||16| . The second says simply that SPPs move from high to low concentration - a 



nonequilibrium "osmotic pressure". These terms are a consequence of the lack of z —>■ —z 
symmetry, and are present in the model of || [7], |Jas well. The third term on the right 
says that in the absence of other forces the director rotates at the same rate as the fluid. 
Together with the fourth term, containing the phenomenological parameter 72, this gives 



rise to the "flow-alignment" [0, [HJ phenomenon, well-known in nematic hydrodynamics. 
These two terms are of course absent in M, |7|, ^] and are central to our work. The ellipsis 
in (|l|) subsumes diffusive terms of higher order in gradients, as well as nonlinear terms. 

The equation of motion for g has the form dtg% = —Vjaij, as in any momentum-conserving 
system, where is the stress tensor. What distinguishes SPP systems is a contribution 
o'fj oc UiHj — (1/3)5^ to the shear stress, as a result of the self-propelling activity. We must 
allow for such a term because no symmetry can rule it out in these driven systems Jll| . To 
derive from the dynamics of individual SPPs note that, by Newton's Third Law, the 
forces exerted by an SPP on the fluid and by the fluid on an SPP are equal and opposite, 



i.e., the combined system of SPPs + fluid is force-free. Thus, in the equation of motion for 
the total momentum density g, the force density associated with each SPP and the ambient 
fluid must integrate to zero. Averaging over the detailed swimming movements, it follows 
that the simplest model for, say, the ath SPP is a rod with axis n a (t), at time t, with point 



forces of equal magnitude / on its ends, directed along ±n Q []nj. For a collection of SPPs 
with centers at r a (t) and ends at r a + ah a and r a — a'h a , where we have allowed for force 
centers asymmetrically disposed about the center of the SPP, this yields a force density 
(divergence of stress) 

-V.*to=fto(r,t) = fJ2^a(t)[5(r-r a (t)-ah a (t)) 

a 

- 6(r - r a (t) + a'h a (t))}. (2) 

The apolar and polar SPPs we mentioned early in this article correspond respectively to 
a = a' and a ^ a'. Expanding the delta-functions in (H) about r a gives 

f (p) „ _ a l±^L fv .Y d h ai h aj 8{x - r Q ) + O(VV) (3) 

In a coarse-grained description this clearly leads to a shear stress 



= ^fc(v,t)(n inj - + O(V) (4) 



as claimed above. Note that only the sum a + a' appears in (Q); to leading order in gradients, 
the self-propelling stress is the same for polar and apolar phases. Using this stress in 
the momentum equation d t gi = —VjUij, Fourier-transforming in space, replacing g by the 
hydrodynamic velocity field u = g/p, and projecting transverse to the wavevector q to 
impose incompressibility (V ■ u = 0) yields, after some algebra, 

^ = -iw q z (\ - 2^±).5n ± - 4" (q±5c) + - (5) 
at q z q z 

where I is the unit tensor, and a ~ fa/ p and w ~ c a are phenomenological constants pro- 
portional to the activity of the SPPs. In (||), nematic elastic torques fl3| and viscous forces, 



both of which are subdominant at small q, have been ignored, as have inertial and other 
nonlinearities. The pressure-like a term despite transverse projection, and the acceleration 
proportional to the bend q z Sn±, arise purely because of self-propulsion |18| . 

Number conservation of SPPs implies d t 5c = — V • j, where the current j = cv oh apart 
from advection by fluid flow and diffusion. To lowest order in gradients and fields, this tells 



us that 

(d t + v d z ) 5c + c v V± ■ 5n ± = (6) 

which says the concentration changes by advection by the mean drift Vq or by local splay 
(since that implies a divergence in the local SPP velocity). 

We can now extract the propagating waves and instabilities mentioned at the start of this 
paper, from the dispersion relation, frequency w as a function of wavevector q, for modes 
varying as exp(zq • r — icjt), implied by equations ([[J), (||) and @. 

Bend-twist waves: Taking the curl of (|l|) and (^j) yields the coupled dynamics of bend 
or twist (V x n_i_) and vorticity (V x uj_). Fourier-transforming in space and time yields 
u = Cbt(0)q, with wavespeeds (see Fig. |l|) 

C6t (0) = ( Cl ±c 2 )cos0 (7) 

where 9 is the angle between the propagation vector and the ordering direction, and c\ and 
C2 (of order the drift speed vq of the SPPs) are phenomenological constants. An analysis 
ignoring hydrodynamic flow, as in || |7|, ||, would instead predict purely diffusive relaxation 
for bend and twist. 

Splay-concentration and drift waves: Taking the divergence of (P and @ results in 
coupled equations of motion for V • ni, V • uj_ , and 5c. The resulting wavelike eigenmodes 
are a generalization of those seen in M. Their speeds obey a messy cubic equation which, 



we can show p0[ , has real solutions for a finite range of parameter values. The speeds as 
functions of direction, for typical parameter values, are illustrated in Fig. 0. 

Viscous effects, and the instability of polar SPP suspensions in the Stokesian 
limit: So far, we have ignored damping terms in the equations of motion since they enter 
at next-to-leading order in wavenumber q, There are three types of such terms: (i) viscous 
damping ~ uq 2 u in the momentum equation (where v is a typical kinematic viscosity of the 
suspension); (ii) Dq 2 5c in the concentration equation, D being a diffusivity; (iii) D n q 2 5n± 
in the director equation where D n is a director diffusivity which we expect to be comparable 
to D. In bacteria, diffusivities are of order 10 -6 cm 2 /s || and can clearly be neglected 
relative to the viscosity. In large-Re systems such as fish, we have no reliable estimate of 
D, but it seems reasonable to assume it is no larger than v. Accordingly, when damping 
terms are included, the modes will remain propagating only for qa < v^ajv = Re, i.e., at 
practically all length scales in the case of fish, but only at exceedingly large length scales for 



bacteria, where Re ~ 10~ 6 to 10~ 4 . The regime of primary interest for bacteria is vq 2 3> voq, 
i.e. qa ^> v^ajv = Re, where viscous damping dominates. For such low Reynolds number 
systems we can use the Stokesian approximation, where the velocity field u is determined 
instantaneously by a balance between viscous and other (in this case self-propulsive) stresses. 
We can thus replace the acceleration in (|5|) by vq 2 w and use it to eliminate u from ([I]) and 
This yields effective equations of motion for the splay V • 5n±_ and 5c, the coupled 
dynamics of which can be seen to yield an unstable mode with growth rate 

r(0) ~ 5(7 2 cos 20+ 1) cos 20 (8) 

with \B\ ~ f<f)/a 2 r] ~ vo<fr/a, <fr being the volume fraction of SPPs. This implies an instability 
just above or just below = tt/4, depending on the sign of B. The mode has a nonzero 
real part ~ ±i7 O <?cos0: this is a convective instability, which is seen if one follows the 
travelling waves. Note that T is independent of the magnitude of q, as a result of the long- 
ranged hydrodynamic interaction in the Stokesian limit. For ga^0 1//2 , director and particle 
diffusion enter to restabilize the mode, so the unstable regime is Re -C qa <C (ft 1 / 2 . 

Instability of SPP nematic suspensions: The dynamics of an SPP nematic suspen- 
sion can be obtained from ([I]), (|5[) and (|G|) simply by omitting all terms not invariant under 
z — > —z, i.e., by setting vq — o\ — a = 0. This can readily be shown to yield the dispersion, 
u 2 ~ —Aq 2 cos 20(72 cos 20 + 1), where \A\ ~ wq . Physically, this can be understood as 
follows. The nematic symmetry z —> —z means that a distortion with wavevector at exactly 
45° has no wavelike response: it doesn't know which way to go, so the 0(q 2 ) contribution to 
the squared frequency w 2 (q) vanishes at this angle, uj 2 must thus change sign as crosses 
45°. For = 71 /A — e, where e is small and has the same sign as A, we thus get an instability 
with u;(q) ~ i\J\e\q — iq 2 . Here the last term appears when diffusive and viscous terms 
and the conventional nematic-elastic torques are included, and will stabilize modes at larger 



wavenumber. For an SPP nematic on a substrate as in [pH| , a wavenumber- independent 
damping enters the momentum equation (|5|) and in general suppresses this instability . 

Fluctuations: In striking contrast to systems at thermal equilibrium, ordered suspen- 
sions of SPPs display giant number fluctuations. To see this, we must include noise, as well as 
the damping terms whose form we have discussed, in the equations of motion. For a general 
SPP suspension one expects noise sources of three kinds: thermal Brownian motion (negli- 
gible for particles larger than a few microns), many-particle hydrodynamic interactions, and 
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intrinsic fluctuations in the self-propelling activity of individual SPPs. In our coarse-grained 
description it is reasonable to assume Gaussian noise sources f M , f c , and f„, delta-correlated 
in time, in the equations of motion for u, c, and n. Momentum and number conservation im- 
ply variances ~ q 2 for f u and f c for wavenumber q — > 0, and nonzero for f„. With these noise 
and damping terms in place we can calculate correlation functions of the various fields in the 
steady state. We find that the structure factor 5(q) = Cq 1 / r exp — zq- t(5c(0)5c(t)) ~ l/q 2 , 
which implies that the variance of the number of particles, scaled by the mean, grows as 
L 2 ~ TV 2 / 3 for a three-dimensional region of linear size L containing N particles on average. 
Physically, this happens because distortions in the director produce mass flow. Since the 
director fluctuations - a broken-symmetry mode - are large, this results in giant number 
fluctuations as well. Similar supernormal number variances were predicted || [7|, || for 
SPPs on a substrate; it is remarkable that they survive here in the presence of long-ranged 
hydrodynamic interactions. 

Let us now estimate the numerical values of the effects we predict. All wavespeeds are 
expected to be of order the drift speed of an SPP (from fim/s to cm/s as one goes from 
bacteria to fish). If we take the kinematic viscosity v to be that of water, then the growth 
rate of the convectively unstable mode in Stokesian polar-ordered SPP suspensions, say for 
bacteria (where velocities are ~ 10/im/s, size a ~ /im, and <fi a few percent), is 0.1s -1 . We 
look forward to experimental tests of our predictions, on living organisms or perhaps on 
artificial SPP suspensions. 

An analysis of the complete (not linearized) equations of motion, now underway, will 
show how nonlinear fluctuation effects renormalize the speeds and dampings of the modes 
presented here, and whether they alter the manner in which number fluctuations diverge 
with size. The full equations of motion will also prove useful in understanding the effect of 
imposed shear flows on the state of order of SPP suspensions, and the role of boundaries 
and topological defects in these systems. These as well as other results on SPP nematics on 
substrates, and an analysis of isotropic SPP suspensions || , will appear elsewhere . 
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